Impact of intrinsic biophysical diversity on the activity of spiking neurons 
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We study the effect of intrinsic heterogeneity on the activity of a population of leaky integrate- 
and-fire neurons. By rescaling the dynamical equation, we derive mathematical relations between 
multiple neuronal parameters and a fluctuating input noise. To this end, common input to het- 
erogeneous neurons is conceived as an identical noise with neuron-specific mean and variance. As 
a consequence, the neuronal output rates can differ considerably, and their relative spike timing 
becomes desynchronized. This theory can quantitatively explain some recent experimental findings. 



In statistical physics, it is often assumed that indi- 
viduals are intrinsically identical. Also in neuroscience, 
identical parameters are typically assumed for all neu- 
rons in the study of neuronal population activity and 
correlation transmission. Real neurons, though, even if 
they are of the same type and located in the same brain 
area, exhibit intrinsic differences. Their morphologies 
and the intracellular concentration of ions, to name just 
two examples, can widely differ, although in principle 
they have been generated by the same mechanisms 
As a consequence, neuronal spike patterns can differ al- 
though neurons receive identical inputs @, H| . Recently, 
in vitro intracellular recordings of isolated mitral cells 
in the mouse olfactory bulb were conducted while they 
responded to identical input [3j (Fig.QJa)). The neu- 
rons displayed diverse output firing rates and pairwise 
correlations. Specifically, the spike correlation coefficient 
obtained with a 1 ms observation window, co- varied with 
the rate difference of the neuron pairs: small differences 
resulted in a wide range of different spike correlations, 
but large differences lead always to small spike correla- 
tion. 

It is of great interest to understand how the biophys- 
ical heterogeneity of a neuronal population contributes 
to neural coding and information processing in neuronal 
networks. Recently, research work has been conducted 
on the coding properties Q and synchronous responses 
5] in a population of heterogeneous neurons. Here we 
investigated more fundamental questions: how neuronal 
heterogeneity can be represented in theory and how it 
can affect the neuronal dynamics and the spiking statis- 
tics in a population of simple leaky integrate-and-fire 
(LIF) neurons using both theoretical analysis and sim- 
ulations. Current models often apply additional inde- 
pendent Gaussian white noises to account for intrinsic 
neuronal heterogeneity. The limitations of this approach 
are addressed first. In contrast, independent noise in- 
put with the same mean and same variance to neurons 
with identical parameters can only capture the trial-to- 
trial variability of neuronal responses, but not the rather 
high diversity in activity and correlation observed in ex- 
periments. Hence, we suggest a more general scheme to 
implement biophysical diversity when either rate or cor- 
relation is of interest. By rescaling the dynamical equa- 
tion, we derive mathematical relations between multiple 
neuronal parameters and the input noise. The main im- 



pact of common input to heterogeneous neurons on rate 
and correlation can be realized by an identical (frozen) 
noise current injection with different values of mean and 
variance, whereas the complete effect is captured by ad- 
ditionally drawing distributed values of the membrane 
time constant and the refractory period. In this scheme, 
the rate difference of heterogeneous LIF neurons can be 
treated analytically. As for correlation, we utilize alter- 
native correlation measures to illustrate that spikes from 
heterogeneous neurons may be desynchronized by sev- 
eral milliseconds, thus escaping the detection by a 1 ms 
observation window. 



I. MODEL 

We consider a population of isolated leaky integrate- 
and-fire (LIF) neurons, each of which has their membrane 
potential V(t) governed by 



r m V(t) = -V(t) + RI(t), 
where the input synaptic current 
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r m = RC is the membrane time constant. R and C are 
the membrane resistance and capacitance, respectively. 
Je (Ji) is the amplitude of an excitatory (inhibitory) 
post-synaptic potential, whereas tj (tk) represents the 
time of the jth (fcth) excitatory (inhibitory) input spike. 
When V(t) — 0, V(t) is reset to V r and a pause for synap- 
tic integration r r is imposed to mimic the refractory pe- 
riod. In the high input regime, the sum of synaptic inputs 
to a neuron can be approximated by a fluctuating input 
noise [10]: 



where 



I(t)=r m (fi + ari(t)), 



[i = Jeve - Jivi, 



° = \IJ% V E 



Jjui. 



(3) 

(4) 
(5) 



rj{t) is a white noise random process such that 
(r)(t)r](t')) = S(t — t'). ve ivi) is the firing rate of the 
excitatory (inhibitory) input. 
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FIG. 1. (a) Correlation coefficient (1ms window) of 589 spike 
trains from mitral cells in vitro receiving an identical input 
as a function of rate difference, adapted with permission from 
1. (b) The same measure for 100 simulated heterogeneous 
LIF neurons using uniformly distributed r m , TV, C and 9 with 
different distribution widths. 

The numerical integration of Eq. Q] in our simulations 
was performed using fourth-order Runge-Kutta method 
with a time step of 0.01 ms. 

II. HETEROGENEITY 

Intrinsic diversity of a population of neurons can be 
directly imposed by drawing neuronal parameters from a 
distribution. Here, we tested the response of 100 isolated 
heterogeneous LIF neurons to an identical fluctuating in- 
put current in the form of Eq. (|3]). Neuronal heterogene- 
ity is implemented by drawing 4 uniformly distributed 
parameters: r m , r r , C (or R) and 9, which, together with 
V r , represent all the independent parameters of a LIF 
neuron in response to a current input. The mean values 
of the uniform distribution are 20 ms, 2 ms, 1 and 1 (in 
arbitrary units) respectively, and V r is fixed to 0. These 
values are used throughout this work. We maintain the 
temporal scale of the dynamics of a typical neuron and 
rescale the potential by setting the mean reset to zero 
and the mean threshold to one. The correlation coeffi- 
cient as a function of the output firing rate difference of 
all possible pairs with different distribution widths (per- 
centage with respect to the mean) from 100 s of simu- 
lations (an example with fi — 0.03 and a = 0.3 shown 
in Fig.[ljb)) highly resemble the experimental findings in 

1 (FigJUa)). 

It is common practice to implement heterogeneity by 
adding independent random inputs to individual neu- 
rons. To test its validity, we provide every identical 
neuron with a common input as the input signal plus 
an independent input with the same statistics among 
neurons to account for heterogeneity, in the form of 
y/c(n + 077(f)) + \/l — c(fj, + a£i(t)) where 77(f) and &(f) 
are independent Gaussian white noises. Fig.[^a) displays 
the raster and the correlation coefficient as a function of 
rate difference when [i = 0.06, a = 0.2 and c = 0.9. The 
rate difference is close to zero and the correlation coeffi- 
cient between any pair is nearly the same Q . Decreasing 
c leads to a drop in spike correlation but has no effect on 
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FIG. 2. Rasterplot, PSTH and correlation coefficient as a 
function of rate difference for the following input into 100 LIF 
neurons: (a) \fc\u, + crr)(t)j + VI — c(/n + cr£j(f)J (c = 0.9); 
(b) fit + <Ti£i(t); (c) our scheme using fa + airi{t) together 
with distributed r m and 7V, consistent with both the experi- 
mental findings 0] and the mathematical analysis. Subscript 
i denotes 'independent'. 



the rate difference. These observations are very distinct 
from both the experimental (Fig.QJa)) and simulation 
(Fig.[IJb)) results. In view of some previous work on 
the reliability of single neurons in response to a repeated 
input 0, [H, [q[ , this implementation may be adopted to 
account for trial-to-trial variability of the same neuron, 
but it fails to account for intrinsic heterogeneity. We 
describe our alternative approach in the next section. 



III. MATHEMATICAL RELATIONS 

In addition to the 5 parameters mentioned above, two 
additional parameters correspond to synaptic inputs: Je 
and J 1 . We analyse the contribution of the heterogeneity 
in the 7 independent neuronal parameters in a population 
of LIF neurons in the high input regime. 

Regarding synaptic inputs, heterogeneity in the pa- 
rameters Je and Jj can be captured by heterogeneity 
in fi and a, according to Eqs. (j4]) and ([5]). For exam- 
ple, if Je is uniformly distributed, /1 is also uniformly 
distributed whereas a 2 is distributed as the square of a 
uniformly distributed variable. 

The other 5 parameters are present in the neuronal 
dynamics irrespective of the type of inputs (current or 
spikes). First, R (or C, depending on the form of writing) 
can be absorbed into /(f) as shown in Eq. (JXJ) so any 
distribution of R can be accounted for by a corresponding 
distribution of /1 and a . 

The difference between 9 and V r , which is the potential 
difference a neuron has to traverse, is a quantity relative 
to the synaptic strengths Je and J/. For instance, lifting 
9, or lowering V r , is equivalent to reducing Je and Jj 
together by the same ratio. Thus, heterogeneity in 9 and 
V r can be included in \i and a by means of rescaling. 

Unlike the above 5 parameters related to the poten- 
tial, the remaining two shaping the neuronal response 
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FIG. 3. Correlation coefficient as a function of rate differ- FIG. 4. Output rate as a function of (a) fx, (b) a and (c) r, 
ence when neurons receive identical noise with uniformly dis- from theory (black) and simulation (colored), 
tributed (a) /i, (b) a and (c) r m separately. 



in the temporal scale, r m and r r , cannot be rescaled or 
captured by /i and a. Their distributions among neu- 
rons have to be accounted for separately. Therefore, in 
the high input regime when the approximation of a fluc- 
tuating input noise is valid, the 7 independent neuronal 
parameters (and their distributions) can be reduced to 
4: jj,, a, r m and r r . Based on this analysis, we sug- 
gest to use distributed values of these 4 parameters to- 
gether with an identical noise 77(f) to account for the ef- 
fects of all the parameters in a population of LIF neu- 
rons receiving identical inputs. This is in contrast to the 
common practice of using independent noises shown in 
Fig.[U[a). We draw the parameters from uniform distri- 
butions fi G [0.015,0.105], cr G [0.1,0.3], r m G [16,24]ms 
and r r G [l,3]ms. In Fig.^c), both the rate difference 
and the correlation coefficient, as well as their relation, 
are consistent with the experimental and our simulation 
results. The significance of this result is tested by replac- 
ing the common noise with an independent noise £i(f) to 
different neuron i. In that case, the correlation coefficient 
lies around zero (Fig.[2jb)). 

The respective contributions of the 4 parameters are 
investigated. Fig. [3] shows the correlation as a function 
of rate difference for 71, a and r m separately. Each re- 
alization is drawn from a uniform distribution of 10 %, 
20 % and 50 % around their mean values, which are 0.06, 
0.2 and 20 ms, respectively. The firing rate of a neuron is 
largely shaped by /j, whereas the distributed values of the 
variance give rise to different degrees of imprecise spiking. 
The wider the distribution of T m and r r , the rate differ- 
ence becomes larger and the correlation drops. When 
r r <C the effect of r r is small (data not shown). 



IV. RATE DIFFERENCE 

The firing rate of a LIF receiving a Gaussian dis- 
tributed noise is known analytically [y, [13] : 



- =T r - T m log(l ), when a = 0; 



(6) 



2 j , 

due u (1 + erf(u)J , otherwise. 



(7) 

Only 6 parameters /i, <r, r m , r r , 9 and V r influence the 
firing rates, of which only the first 4 are independent. 



Changing any of them can result in a rate difference as 
shown in Fig. 31 and this explains the wide distribution 
of firing rates in a heterogeneous neuronal population. 



V. IMPRECISE SPIKING 

The raster plot in Fig.^c) shows population synchrony 
with spike-time jitter. Low average correlation coeffi- 
cients in a population of neurons do not necessarily im- 
ply an asynchronous state of the population. Removing 
spikes from one of the two identical spike trains can also 
reduce the correlation coefficient significantly. When we 
compute the correlation coefficient in our data with a 
larger bin size, larger values for the the correlation co- 
efficient are obtained, as shown in Fig.[SJa). This is be- 
cause some spikes become 'coincident' only for larger bin 
sizes. The output spike trains behave like jittered spikes 
as discussed in earlier theoretical studies [11] . Whether 
the decorrelation due to neuronal heterogeneity is signifi- 
cant depends critically on the bin size, or the integration 
window of the neurons receiving such inputs. 

In view of a significant number of spikes jittering out- 
side a 1 ms bin, we look into the cross-correlation function 



r X y(r) 



c xy (r) (x(t)y(t + T ))-(x(t))(y(t)) 



(8) 



where x(t) and 77(f) denote two output spike trains in dis- 
crete time, consisting of and 1 with bin size of 0.1ms. 
x(t) is assigned to be the spike train with lower spike 
count. c xv (t) is the covariance function and cr x , a y de- 
note the standard deviation of the two spike trains, con- 
sidered as discrete signals. Fig.JSJb) shows that the mean 
r xy (r) over all pairs is positive in the small neighbour- 
hood of t, indicating a higher than chance level to observe 
spikes. Spikes arc jittered, instead of asynchronous. In 
addition, r xy {r) is asymmetric, and its area is skewed to- 
wards negative r, indicating that the higher firing rate 
neuron is more likely to lead in terms of spiking [l2|, [Hj] • 
We further look at the normalized total cross- 
covariance k [ill, E3 



\] I-00 drc xx {T) dT'c V y(T f ) 



(9) 



which is an overall measure for the fraction of the spikes 
that are correlated above chance level. Fig.[Stc) shows 



4 




Rate difference (Hz) 



FIG. 5. (a) Correlation coefficient as a function of window 
size (colored lines indicate 5 examples; the black line indicates 
the average over all pairs), (b) Cross-correlation. Positive 
correlation at negative r indicates that the higher-rate neuron 
is leading, (c) k, which takes jittered spikes into account, as 
a function of rate difference. 



that k is quite close to unity, and show a weak depen- 
dence on rate difference. This indicates that spikes would 
not be decorrelated when a larger time window is consid- 
ered. 



VI. DISCUSSION 

We remark that the plot of correlation coefficient as a 
function of rate difference from our simulations (such as 
Figs. Hlb) and[2fc) can fit the experimental results in [3[ 
more satisfactorily by introducing a random delay of up 
to 1 ms to every incoming input spike. The spike correla- 
tion distribution then becomes broader, and the regime 
of low correlated output at small rate differences can also 
be reached (data not shown) . Such effect could be due to 
some (unknown) biological mechanism not captured in a 
simple LIF neuron model. 

We emphasize that common input into heterogeneous 
neurons is better realized by a shared noise with dis- 
tributed mean and variance, than by the common prac- 
tice of adding independent noises of the same statistics. 
This insight is based on both the mathematical analy- 



sis presented here and the in vitro experimental findings 
Q. As far as firing rates and spike correlations are con- 
cerned, a distribution of mean and variance of the input 
account for most of the experimental observations. Di- 
versity of r m and r r have a smaller effect of the quantities 
in question, nevertheless, including them can account for 
the full degree of heterogeneity. 

In the raster plot of a neuronal population with sim- 
ilar rate differences and spike correlations, synchrony is 
obvious, although spike times are not precise. Spikes, if 
present, are jittered in the millisecond range, which can- 
not be captured by the 1 ms temporal window used for 
analysis. This is why using larger bins leads to larger 
values for spike correlation. We emphasize that neu- 
ronal heterogeneity alone does impose an appreciable 
decorrelation effect on the population activity. However, 
whether decorrelation is functionally significant depends 
on the readout of the downstream neurons. On top of 
that, a network of heterogeneous neurons may give rise 
to richer network dynamics. It remains to be explored 
whether such intrinsic heterogeneity can facilitate other 
decorrelation mechanisms to increase the amount of in- 
formation flow [1, [TH, [l6[ . Given the significant reduc- 
tion in spike correlation among heterogeneous neurons, 
research concerning correlation transmission must take 
neuronal heterogeneities into consideration. 
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